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^ ■ Abstract 
in 

We calculate the third coefficient of the lattice (5 function in pure Yang- Mills 
'■ theory. We make use of a computer code for solving perturbation theory 



analytically on the lattice. We compute the divergent integrals by using a 
method based on a Taylor expansion of the integrand in powers of the external 
Jjn | momenta in 4 — e dimensions. Our results are in agreement with a previous 



X 



calculation by M. Liischer and P. Weisz where the authors used a different 



technique. We also show how this new coefficient modifies the scaling function 
on the lattice in both the standard and energy schemes. In particular we show 
that asymptotic scaling is extremely well achieved in the energy scheme. 
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I. INTRODUCTION 



The relationship between the bare coupling and the cutoff is an essential ingredient in 
lattice calculations. Usual Monte Carlo simulations require the knowledge of this relationship 
far from the critical point, where corrections to asymptotic scaling could become relevant. 

In order to check the relevance of these scaling corrections, in this paper we present a 
calculation of the first non-universal coefficient of the lattice (3 function in pure Yang-Mills 
theory. This renormalization group function can be written as 

P L (go) = -o^ U= -6fog - %gl - fal - • • • (i-i) 

where a is the lattice spacing, g the lattice bare coupling, g r is the renormalized coupling 
constant and \i the subtraction point. It establishes how the bare coupling and the cutoff a 
must simultaneously vary to keep the renormalized quantities fixed. The first two coefficients 
in eq. (1.1) are known as they are scheme-independent. Our purpose is to compute the third 
coefficient, b\ . 

In a recent paper the authors calculated a quantity related to this coefficient using 
a coordinate space method [[J for evaluating the lattice integrals. In our calculation we 
have used a different technique where firstly the integrands are Taylor expanded in powers 
of the external momenta and then computed term by term with the introduction of an IR 
regulator. The idea is a generalization to two loops of the procedure introduced in |J. Our 
result is an independent check of the calculation in reference |IJ. O ur integration method 
can be straightforwardly generalized to the case with fermions. 

In order to facilitate comparison we have adopted the notation of ref. [[I]]. 

The plan of the paper is as follows. In section 2 we will show the method we followed 
to determine b\. We will give explicit formulae that relate the coefficient b\ to the same 
coefficient in the renormalized f3 function (see eq. (2.9)). In section 3 we give the final result 
of our calculation as the two-loop contribution to the coupling renormalization constant 
and we derive the value for b\ and discuss its consequences for the scaling behaviour of 
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dimensionful quantities in present-day Monte Carlo simulations for both the standard and 
energy schemes. Finally, in section 4 we introduce our integration method by solving a 
typical two-loop integral. The conclusions are presented in section 5. A list of superficially 
divergent integrals used in our computation are shown in the Appendix. 

The involved algebra of the lattice perturbation theory was carried out by making use 
of a computer code. The main lines of this code were explained in [|J and used to compute 
the three-loop perturbative background of the topological susceptibility in pure Yang-Mills 
on the lattice and the three-loop lattice free energy in the same theory, || . For the purpose 
of the present paper, this code was extended to include form factors. 



Although &2 is a three-loop quantity, it can be computed by evaluating two-loop diagrams 
on the lattice provided the corresponding coefficient in the renormalized f3 function, 62, is 
known. The renormalized f3 function is defined as 



This function depends on the renormalization scheme only, hence it must be the same for 
any regularization. In this section we will show how the coefficient h 2 can be derived from 
this f3 function. 

The bare coupling constant on the lattice go and the renormalized one g r (for instance 
in the MS scheme in the continuum) are related by 



II. A TWO-LOOP CALCULATION 





g = Z(g ,fia)g r 



(2.2) 



where 



Z (g , fia) 



1 + Z l9 l + Z 2 g 4 +■■■ 



Zi = Zio + Zu In /id 



Z 2 = Z 2 o + Z 2 i In \xa + Z 2 2 In 2 fxa 



(2.3) 
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The f3 L function (1.1) can be computed from the previous renormalization constant by 
calculating the derivative of (2.2) with respect to the lattice spacing a 

<*(™^fiW (2.4) 
Z(g ,jja) \ dg da J 

This equation yields 6q = Z n , b\ = Z 2 \ + Z w Zu; while the fmiteness of the beta function 
implies Z\ x + 2Z 22 = 0. 

Now, combining the derivative of eq. (2.2) with respect to \i 

= (3(g r )Z(g , H + g ^ 9 ^ (2.5) 

with eq. (2.4) we obtain 

oL ( \ P(9r(g ))Z{g Q ,iia) 

P i90) = - dlnZ(g ,,a) - ^ 
d In g 

This equation relates the perturbative coefficients in eqs. (2.1) and (1.1). In particular we 
obtain 

60 = bo, (2.7) 

61 = h, (2.8) 
b L 2 = b 2 - 2hZ 10 + b Z 2 10 + 2b Z 20 . (2.9) 

Eqs. (2.7) and (2.8) establish the well-known result that the first two coefficients are inde- 
pendent of the renormalization scheme. Eq. (2.9) implies that if we know the value of b 2 in 
some renormalization scheme, we need only to perform a two-loop calculation in that scheme 
on the lattice to obtain b 2 . Throughout this work we will use the MS scheme. Had we used 
eq. (2.4) then a lattice three-loop calculation would have been needed. The coefficient Z w 
was computed in reference 0, its value for the gauge group SU(N) being 

z » = N iih + m;-h-fr-T*Y (2 - 10) 

In this expression Pi and P 2 are numerical constants, Pi = 0.15493339, P 2 = 0.024013181. 
The coefficients bo and b\ are 



11 N 

bo 



3 16vr 2 ' 

2 



34 / iV . 

bi = ^T T7-,) • 2.11 



3 V16tt 2 

Moreover, the coefficient 62 i n the MS scheme is || 

2857 / N \ 3 . , 

62 "(1^) ■ (2 - 12) 

The computation of the coupling renormalization constant, eq. (2.3), is easier in the 
background field gauge. In fact, this renormalization constant has a simple relationship 
with the background field renormalization constant 

Z A (g , fia)Z(g , fia) 2 = 1 A^ = Z A (g , ^a) l ' 2 A r ^. (2.13) 

Relation (2.13) has been recently proven also on the lattice ||. Hence we need to compute 
only a two-loop self-energy of the background field on the lattice. In our computation we 
used the Wilson action in the background field gauge ||. In this formulation the links are 
written as 

U^x) = U«{x)U%{x), U*(x) = e i9oQ » {x \ U c J(x) = e ia9oA » {x \ (2.14) 

where Qn(x) = T c Q c ^(x) and A^{x) = T c A c tl {x) are the quantum and background fields 
respectively. The following gauge-fixing term preserves gauge invariance of the background 
field 

D-Q v (x) = U*~ x (x - jj)Q v {x - p)U*(x - A) - Qu(x). (2.16) 

We worked in the Feynman gauge, a r = 1. This gauge fixing produces the following Fadeev 
Popov action for the ghosts fields u and u; 



-^9o [Q^), [Qn(x),D+u(xj\]+---), (2.17) 
D^(x) = U c J(x)uj(x + p)U$-\x) - u(x). (2.18) 
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Finally the change of integration variables from links to vector fields yields a jacobian that 
can be rewritten as the measure action 

S m = ±-Ngl E E Tr QMQM + ■■■ (2-19) 

X fl 

In all expansions, eqs. (2.17) and (2.19), we have written only the relevant terms for our 
two-loop computation. 

The full action is therefore S = Swuson + S g f + S g h + S m . The vertices needed in 
perturbation theory are obtained as usual by expanding the exponential of the previous 
action in powers of the fields. In ref. 0] we described in some detail our automatic procedure 
for generating vertices. 

The Feynman diagrams contributing to the background gluon self-energy are shown in 
Figure 1. The bare and renormalized self-energies of the background field are related by eq. 

(E3) 



(1 - u r (p, fi, g r )) = Z A {1 - u(p, a, g )). (2.20) 

In the computation we used the bare quantum fields. Therefore the gauge parameter must 
be explicitly renormalized. Diagrams 1 — 4 correspond to this renormalization. Up to two 
loops the bare self-energy can be written as 

v(p,a,g ) = glv {l) + g^ 2) + a r (z% + Z^lnfia) ^-\ ao=olr , (2.21) 

where the relation cto = Z®a r has been used. Z® is the quantum field renormalization 
constant. We need this constant at one-loop, Z Q = 1 + g^iZ 1 ^ + Z^ In fia) 0, 

Zll ""YToV' 

Z% = N ( —P x + - l — l — + -P 2 ) . (2.22) 

10 V72 1 16 8N 2 48tt 2 3 J K J 



The last term in eq. ( |2.21| ) corresponds to the diagrams 1 — 4 in Figure 1. 



The renormalized self-energy can be written similarly as u r (p,fi,g r ) = gfor + 
g^v^ + ■ ■ ■■ By using this relation and eqs. (2.3), (|2.13|) , (|2.21|) we obtain 



2Z 2 + (Z,) 2 = - - (Z% + Zg In //a) a r ^. (2.23) 

In this expression, we took advantage of the gauge parameter independence of the coupling 
renormalization constant Z in the MS scheme to substitute by in the derivative [10 



The explicit expressions for the one-loop and two-loop renormalized self-energy of the 
background field were computed in ref [11]. In the MS scheme they are 

m iV / 11 p 2 205 3 1 „\ , 

^ 2, = (i^) 2 (- 81n ^ + ^- 6 « 3) )' (2 ' 25) 

Here £(3) = 1.202057 is the Riemann's zeta function and fi is the renormalization point. 
The two-loop self-energy (2.25) is written in the Feynman gauge. 



III. THE RESULTS 



The result for each diagram can be written as A p {p)A a {— p)T pa (p) with T pa (p) some 
tensor. Leaving out the fields and summing over p, a we obtain TrT(p). The gauge invariance 
of the final result guarantees that the sum of all diagrams yields 3p 2 (1 — u(p, a, go))- 

For some diagrams we have computed the whole expression A p (p)A a (—p)T pa (p) for pur- 
poses of checking. In all cases we obtain agreement (within the precission used) with [0. 
Therefore we do not show the results for each diagram. 

After collecting all individual contributions we get for the following expression 

v ® = i n ( a y ) + + Cl + c 2 N 2 

32tt 4 v v ' 128 iV 2 

ci = -0.01654 

c 2 = 0.007230. (3.1) 

Now, using eq. (2.23) we get the two-loop term of the coupling renormalization constant 
Z(go,/j,a) in eq. (2.3) 
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1 \ 2 / 340.9 „ 
■ 358.5 — 188.7 Af + 



16tt 2 / V ' N 2 
N 2 



/193 11 2 110 ,„ 968 2 \ 11 2 \, 121 iV 2 , 2 \ , n 

+ — 7T 2 + — 7T 2 P! + 7T 2 P 2 7T 2 111 UO 111 2 UO . 3.2 

V 18 6 27 9 V 3 J p 18 p J K ' 

Notice that as predicted Z\ x + 2Z 22 vanishes. The coefficient of the logarithm in the previous 
equation contributes to the second coefficient of the (3 function, b x . The other term in eq. 



[2D together with eqs. (2.9) and (2.10-2.12) gives the desired result for b 
,r / N \ 3 / 22 vr 4 1 / 2 2816 4 . 

&2= TTT^ -^7^-^ 87T 2 + — 7T 4 Cl + 



L 

2 ! 



2 U6ttV V ^ 4 N 2 

1522 2 2816 , 124 „ 1408 9 w A . , 

— + 4 vr 2 - — vr 4 c 2 + — vr 2 P 1 + — tt 2 P 2 - 22 C (3) (3.3) 

and numerically 

, T ( N \ 3 / 1433.8 2143. \ , , 

Now, by integrating the lattice j3 function, eq. (1.1), we get the following dependence of the 
lattice cutoff on the bare coupling 

aA L = exp (-1/2M 2 ) (b g 2 )- b ^ 2b o (l + ± (b\ - b%) g 2 + • • \ , (3.5) 

where is the renormalization group invariant mass in the lattice regularization. For the 
gauge group SU (3) eq. ( |3.5| ) becomes 

/I6vr 2 \ 51/121 / 8tt 2 \, 2 \ 

aAL = \Tigi) exp yWo ) ( + °- 189% °) • (3 - 6) 

Thus in the region go w 1 the last factor in eq. (|3.6| ) brings about a ~ 20% correction. In 



the SU(2) gauge group this correction amounts to 0.08324 g 2 .. At large A" the correction 
behaves linearly as 0.1048 N g 2 ,. 



We can also calculate the correction to asymptotic scaling in the scheme of the energy |T2 
The weak expansion of the average plaquette is [|5|,|l~3" 

( 1 - ^ D ) = ™i£o + w 2 gt + w z gl + ■■■ 

N 2 - 1 N 2 -I f 1 
wx= 77- w 2 = 0.02043 - 



! A" 4 V ' 32 A" 2 

N (N 2 - 1) / 0.01812 0.01852X , . 

<*= — 0.006354 - + . (3.7) 



The coupling constant in this scheme can be defined as 

With this definition, the third coefficient of the beta function is 



and the scaling function in SU (3) becomes 



-, . , 2 \ 51/121 / 2 

167T \ / 87T 



aA^^j exp^j (1 + 0.01161,2 
A B = exp (g^) A ^ = 2 - 0758 ( 3 - 10 ) 

In SU{2) the correction is 0.003530 g%. 

The previous results show that as expected the first correction to asymptotic scaling in 
the scheme of the energy is much smaller than in the usual standard scheme. This fact 
indicates that in the scheme defined by eq. fl3.8|) asymptotic scaling is possibly satisfied 
within a few per cent at bare couplings as large as g% ~ 1. 



IV. THE INTEGRATION METHOD 



In this section we introduce our integration method. First we show how we computed 
UV-divergent integrals and then we will say a few words about the convergent integrals. 
The method for UV-divergent integrals will be shown by solving the following example 



71 d A k r +7T d 4 q 1 

ToZyI / /r>_N4 -2-2 - 2 2 

tt [ATT) J-tt [ATI) q kq + kk + ap 



E (ap) = I 77^X1 I ,0X^2 ^,2, — 2 - (4-1) 



Here p is the external momentum and q = 4X^sin (q fJ- /2). This integral diverges when 
a — > 0. Similar one- loop divergent integrals can be computed by the following technique 
(see also ref. Jl4|]): ^ define the same integral in D > 4 dimensions; Taylor expand 
this integral in the external momenta until the remainder of the expansion is UV-finite; in) 
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integrate all terms in this expansion including the remainder and taking into account the 
IR-divergences that are generated. The remainder can be computed in the continuum as it 
is UV-finite. Eventually all IR-poles must cancel giving an IR-finite expression. Notice that 
the computation of the remainder in the continuum requires an exchange of the limits a — > 
and em — > 0; this delicate operation is valid provided adequate subtractions are made on 
the integrand. 

The method we used to compute the integral (4.1) is a generalization of this idea beyond 
one loop. We must simply apply the same steps to every loop. Then the integral in eq. (4.1) 
becomes (k is an arbitrary mass scale) 



w x , M f f +7T d D k r+* d D q ( 1 



I 1 1 ^ 



+ k (<f 



yk k + ap \k 2 J 



+ (4-2) 
\ 



, * d D k r+* d D q I 1 1 



-7T (2n) D J-n (2n) D \kk + ap{q 2 f qq + k(k 2 ) 2 (k 2 )\q 2 f j 
where we put D = 4 — 2e (from now on we will omit the subscript IR to the e). Now the 
integral in (4.2) is UV-finite and we can exchange the UV and IR limits to rewrite it as 
f d^_ r d D q p(p + 2k) k(k + 2q) 

V (2tt) d Jr° (2tt) d (k + p) 2 (k 2 ) 2 (k + q) 2 (q 2 ) 2 { h { ' 

which can be computed by continuum methods. 

The first and third terms in eq. (4.3) can be combined to give 

4e f+* d D q 1 f+« d D k ( 1 1 \ 

{Ka) L (2n)°(q 2 ) 2 L [2,r \k 2 kTap 2 (k 2 ) 2 ) " ^ 

The /c-integration is UV-finite, therefore we can solve this integral by perfoming the corre- 
sponding continuum limit. 

The second term in eq. (4.3) is IR-divergent and can be calculated by performing a chain 
of adequate additions and subtractions. We can rewrite it as 

/ ^A f +7T d D k 1 /■+* d D q ( 1 1 \ . . 

( "°> L j^WrJ- l^T? " Wr) + (4 ' 6) 

[+<■ d D k r+» ,\ D q 1 1 
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This operation has not modified the bad IR and UV properties of the integrand. However, 
eq. ( |4.7| ) can be computed numerically as we will show below. Let us study eq. (|4.6|) . We 
would know how to solve this integration if it was written in the continuum. Following this 
hint we add and subtract the analogous expression in the continuum for the integration in 
q to get 

4e /■+*■ d D k 1 d D q I 1 1 1 1 \ 

4e /•+* d D k 1 /•+*■ d D q I 1 1 \ 



As expected, eq. (|4.8|) is IR-finite so that it can be computed numerically in 4 dimensions. 
Furthermore, the subtraction of fractions in the g-integration of eq. ( f4.9|) makes its k- 
integration IR-finite. Therefore all IR-divergences in eq. (|4.9| ) are due to the q 2 ~ region 
in the g-integration. Hence a possible way to solve this integral is to split the g-integration 
into two regions 

d D k 1 r d D q ( 1 1 \ 

M J — 7r (2tt) d (jfe2)2 J r d (2^p \ q 2( q + k) 2 ~~ ~~ ( ' ' 

4e r-Hr d D fc 1 /• d D g / 1 1 \ 

where IR D — ii D means that the integration domain covers all the euclidean space-time except 



for a hypercube of side ir around the origin. Now, eq. ( 4.11]) is IR-finite and can be evaluated 



numerically in D = 4 dimensions. Integral Q4.10 ) can be solved by continuum techniques 



Some one-loop IR-divergent integrals are needed to 0(e). They are 



W^d = TT-2 I" + 1-838199 + 2.283043 e + 0(e 2 ' 
J-TT (27r) w (q 2 ) z 16n 2 \ e 



C (& W? (q2) " = Tt?(-i + 3315103 " 36121 f + 0(£2) ) ■ (412) 

The first one can be evaluated by using Schwinger parameters. The second one needs first 
an expansion in powers of e, then the Schwinger parametrization can be used. 

After evaluating the integrals in eqs. ( |4 . 4|J4 . 5| , P~T| , |4 . 8| , |4TTD| J4 . 1 1| ) and adding up the partial 
results we obtain 
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E ( a P) = (r^) (28.01 - 6.7920 1n(a 2 p 2 ) + iln 2 (a 2 p 2 )). (4.13) 

As expected, this result is independent of the IR regularization we used. 

In table I we show a list of the main two-loop divergent integrals used in this work. All 
entries in this table have been computed by using the above technique and the results agree 
with ref. ||. 

All the other superficially divergent integrals appearing in our computation can be rewrit- 
ten in terms of this basic set plus some finite integrals. The finite lattice integrals have been 
solved either by using a method based on Schwinger parameters explained in ref. p],[15|j or 
by an extrapolation to infinite lattices from their values on finite lattices. The extrapolating 
function was chosen to be r + J2nm r nm{^ m L)/L n where L is the lattice size and n and m 
free parameters. The sum on n and m runs over a few integer values (typically m = — 1 
and n = 2). The infinite lattice value of the integral is r . Both methods agree. We did 
not push the precision to more than four or five digits (for practical purposes in numerical 
simulations, knowledge of four digits is sufficient). 

In the appendix and in table II we give the list of supercifially divergent integrals en- 
countered during our computation. 

V. CONCLUSIONS 

We have calculated the third coefficient of the lattice (3 function for the Yang-Mills theory 
with gauge group SU(N). A computer code was used to simplify the algebra of perturbation 
theory on the lattice. We have also introduced a method to compute the two-loop lattice 
integrals. The method can be generalized to any number of loops and with the inclusion of 
fermions. Our result agrees with a similar computation of reference Q. 

The corresponding corrections induced by this new coefficient cannot be neglected near 
the scaling window for the standard scheme. For instance the modification in the 577(3) 
gauge group amounts to ~ 20% when g w 1. Instead, in the energy scheme defined by eq. 
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( P-8| ) the correction to asymptotic scaling in SU(3) is of the order of ~ 1% at qe ~ 1- The 
improvement in SU(2) is also significant: a ~ 8% correction to asymptotic scaling in the 
standard scheme versus ~ 0.3% in the energy scheme when gg « 1. This result emphasizes 
the role of the energy schemes in numerical simulations which require good asymptotic 
scaling. 



VI. APPENDIX 

We present in table II a set of superficially divergent two-loop integrals encountered in 
our computation. They can be written in terms of the basic set in table I plus finite integrals 
by using algebraic manipulations based on the identity 

rTs 2 = r 2 + s 2 +2^2 M^) sin^) - 8 £ sin 2 (^j sin 2 , (6.1) 

valid for any internal or external momenta r, s. 

This table together with table I is a basis for all superficially divergent integrals with 
gluonic propagators, depending on one external momentum. We evaluate these integrals 
in the continuum limit. Integrals with more powers of momenta can be evaluated either 
by Taylor expansion or, in the presence of sub divergences, by subtractions as in section 4, 
leading to a Taylor expandable part plus a product of one-loop integrals. 

Some trivial cases have been omitted (cases with an inverse propagator in the numerator, 
and cases containing k^C^ v {k) = 0, see M). 

The notation is as in table I: q, k are internal momenta and p is external. On the other 
hand, (r • s) = srn ?V srn s fi an d 

-w d 4 k /■+*■ d 4 g (r ■ s) 



F rs (ap) = J —j_ 



'-7T (27r) 4 J-w (27r) 4 g2£ + (? Ip-k^p 

j . r +7T d 4 k r+ 7T d 4 q (r ■ s) (u- v) 

&rsuv\ a 'P) = / 77m / /r> \A 2 - - — - 2 

(2tt) 4 J-* (2tt) 4 q 2 k+q [k 2 fk+p 



w d 4 k r +7T d 4 q (r • s) (u • v) 

tt (27r) 4 J-tt (27r) 4 q 2^ q 2 p,j^p ^Zp 2 



a,„M-/ tHi/ t^i -^J-s-, - (6-2) 
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Figure captions 

Figure 1. Set of Feynman diagrams contributing to the two-loop self-energy of the background 
field gauge. Straight, dashed and wavy lines correspond to quantum fields, ghosts and 
background external fields respectively. The black circle stands for the measure vertex 
and the black square for the counterterms. 



Table captions 

Table I. Basic two-loop divergent lattice integrals. All expressions in the second column are in- 
tegrated with the measure d 4 q/(2n) 4 d 4 A;/(27r) 4 . The constant C(3) = 1.202057 
is the Riemann's zeta function. 

Table II. Superficially divergent integrals. The notation follows the appendix. Pi and P 2 have 
been defined in the text. 
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Table I 



integrand 


result 


1 

-2-2 -— ~ 2 
q k q+k+ap 


(^) 2 (100.8 + a 2 p 2 (-1.857 + I ln(aV))) 


1 

-2-2- -2 - — - 2 

q k q+k k+ap 


(J 2 ) 2 (28.01 6.7920 ln(aV) + J ln 2 (aV) 


) 


1 

-2-2— — •! a 

g fc g+fc q—ap k+ap 


(l£*) aV 6 C(3) 


sinfejj sinfc„ 

-2 /-2\ 2 2 2 

g 1 k 1 g+fc k+ap 


(l^) 2 (^ (2-646 -|ln(aV)) + 
(7.27 '4 - 3.146 ln(a 2 p 2 ) + ± ln 2 (a 2 p 2 )) ) 




sing M sing„ cos(fc M /2) cos(fc„/2) Ny 


(irfftT (0-07838 -£ln(aV)) + 


) 


/-2\ 2 2 X 

1 A; 1 k+ap 

{(q + ltq)- 1 ~ (g 2 ) j 


|V (-1.259 + 0.4435 ln(aV) - ^ln 2 (a 2 p 2 )) 



Table II 



integral 


result 


Fpk(ap) 


(lebO 2 («V (" 6 - 774 +(! + 8vr 2 P 2 ) In(aV) - ^^)) 


F pq (ap) 


(i^) 2 («V (1-613 + (-§ + £ (Pi " 8P2)) In(aV) + "S^)) 


Fqk{ap) 


(^) 2 (-9.246 + aV (0.9368 - ^)) 


Eppkq(ap) 


(i6^) 2 («V (- 4 -286 + (| - vr 2 (P x - 8P 2 )) ln(a 2 p 2 ) - >^>)) 


E p kpq(ap) 


(ife)' (*V (-1-490 + (A _ ( Pl _ 8 P 2 )) In(aV) - ^)) 


Epkkq(ap) 


(^) 2 (a 2 p 2 (1.097 +(-§ + £ (Pi - 8P 2 )) ln(aV) + ^1)) 


G PP pk{ap) 


(^) 2 (-3C(3)a 2 p 2 ) 


Gp pk q(ap) 


(le^)Vp 2 (-4-214 + ln(aV))) 


Gpk P k( a P) 


(iB^) 2 («V (5-337 +(-! + £ (Pi - 8P 2 )) ln(aV) + ^V))) 


Gpk pq {ap) 


(^) 2 (a 2 , 2 (-2.732 + ^)) 


Gpkkq(ap) 


(^) 2 (aV (2.846 + (-g + £ (Pi - 8P 2 )) In(aV) + ^)) 


G k qkq(ap) 


(^J 2 (6.486 + aV (0-1405 ^V))) 
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Figure 1 




Figure 1 (cont.) 



